to remove/jagsmodelTOrun.R

mydoseresNMA <- function(){
  # Begin Model Code

  s.beta.2[ref.index] <- 0

  s.beta.1[ref.index] <- 0

  for(i in 1:ns){ # Run through all NS trials
    rr[i,1]~dbinom(p0[i],nn[i,1])
    delta[i,1] <- 0

    #DR[i,1] <- 0 # Dose-response model is 0 for baseline arms
    mu[i] ~ dnorm(0,0.001)

    for (k in 1:na[i]){ # Run through all arms within a study

      rhat[i,k] <- p[i,k] * n[i,k]
      resdev[i,k] <- 2 * (r[i,k] * (log(r[i,k]) - log(rhat[i,k])) + (n[i,k] - r[i,k]) * (log(n[i,k] - r[i,k]) - log(n[i,k] - rhat[i,k])))
      r[i,k] ~ dbin(p[i,k], n[i,k])
      logit(p[i,k]) <- theta[i,k]
      theta[i,k] <- mu[i] + delta[i,k]
    }

    resstudydev[i] <- sum(resdev[i, 1:na[i]])

    for(k in 2:na[i]){ # Treatment effects

      delta[i,k] <- DR[i,k]

      DR[i,k] <- ((s.beta.1[t[i,k]] * dose1[i,k]) + (s.beta.2[t[i,k]] * dose2[i,k])) - ((s.beta.1[t[i,1]] * dose1[i,1]) + (s.beta.2[t[i,1]] * dose2[i,1]))
    }
  }

  for (k in c(1:(ref.index - 1), (ref.index + 1):ndrugs)){ # Priors on relative treatment effects
    mult[1:2,k] ~ dmnorm(d.prior[], inv.R[1:2, 1:2])
    d.2[k] <- mult[2,k]

    s.beta.2[k] <- d.2[k]

    d.1[k] <- mult[1,k]

    s.beta.1[k] <- d.1[k]
  }

  totresdev <- sum(resstudydev[])
  Omega[1,1] <- 1
  Omega[2,2] <- 1

  for (r in 1:2) {
    d.prior[r] <- 0
  }

  inv.R ~ dwish(Omega[,], 2)

  for (r in 1:(2-1)) {  # Covariance matrix upper/lower triangles
    for (c in (r+1):2) {
      Omega[r,c] <- 0   # Lower triangle
      Omega[c,r] <- 0   # Upper triangle
    }
  }
  for (m in 1:ns){ ## for each study
    #   # rr[m,1] ~ dbinom(p0[m],nn[m,1])
    logit(p0[m]) <- OO[m]
    OO[m] ~ dnorm(Z, prec.z)
    #   # z[m] ~ dnorm(0,0.001)
  }
  # priors
  Z ~ dnorm(0, 0.001)
  prec.z <- pow(sigma.z,-2)
  sigma.z ~ dunif(0,10)

  for (k in c(1:(ref.index-1),(ref.index+1):ndrugs)){
    for( j in 1:nd.new[k]){
      OR[j,k] <- exp(s.beta.1[k]*new.dose[k,j]+  s.beta.2[k]*f.new.dose[k,j])
      odds.drug[j,k] <- OR[j,k]*exp(Z)
      p.drug[j,k] <- odds.drug[j,k]/(1+odds.drug[j,k])
    }
  }
}



mydoseresNMR <- function(){
  # Begin Model Code

  s.beta.2[ref.index] <- 0

  s.beta.1[ref.index] <- 0

  for(i in 1:ns){ # Run through all NS trials
    rr[i,1]~dbinom(p0[i],nn[i,1])
    delta[i,1] <- 0

    #DR[i,1] <- 0 # Dose-response model is 0 for baseline arms
    mu[i] ~ dnorm(0,0.001)

    for (k in 1:na[i]){ # Run through all arms within a study

      rhat[i,k] <- p[i,k] * n[i,k]
      resdev[i,k] <- 2 * (r[i,k] * (log(r[i,k]) - log(rhat[i,k])) + (n[i,k] - r[i,k]) * (log(n[i,k] - r[i,k]) - log(n[i,k] - rhat[i,k])))
      r[i,k] ~ dbin(p[i,k], n[i,k])
      logit(p[i,k]) <- theta[i,k]
      theta[i,k] <- mu[i] + delta[i,k] +(gamma[t[i,k]] - gamma[t[i,1]])*cov[i]
    }

    resstudydev[i] <- sum(resdev[i, 1:na[i]])

    for(k in 2:na[i]){ # Treatment effects

      delta[i,k] <- DR[i,k]

      DR[i,k] <- ((s.beta.1[t[i,k]] * dose1[i,k]) + (s.beta.2[t[i,k]] * dose2[i,k])) - ((s.beta.1[t[i,1]] * dose1[i,1]) + (s.beta.2[t[i,1]] * dose2[i,1]))
    }
  }

  for (k in c(1:(ref.index - 1), (ref.index + 1):ndrugs)){ # Priors on relative treatment effects
    mult[1:2,k] ~ dmnorm(d.prior[], inv.R[1:2, 1:2])
    d.2[k] <- mult[2,k]

    s.beta.2[k] <- d.2[k]

    d.1[k] <- mult[1,k]

    s.beta.1[k] <- d.1[k]
  }

  totresdev <- sum(resstudydev[])
  Omega[1,1] <- 1
  Omega[2,2] <- 1

  for (r in 1:2) {
    d.prior[r] <- 0
  }

  inv.R ~ dwish(Omega[,], 2)

  for (r in 1:(2-1)) {  # Covariance matrix upper/lower triangles
    for (c in (r+1):2) {
      Omega[r,c] <- 0   # Lower triangle
      Omega[c,r] <- 0   # Upper triangle
    }
  }
  for (m in 1:ns){ ## for each study
    #   # rr[m,1] ~ dbinom(p0[m],nn[m,1])
    logit(p0[m]) <- OO[m]
    OO[m] ~ dnorm(Z, prec.z)
    #   # z[m] ~ dnorm(0,0.001)
  }
  # priors
  Z ~ dnorm(0, 0.001)
  prec.z <- pow(sigma.z,-2)
  sigma.z ~ dunif(0,10)

  for (k in c(1:(ref.index-1),(ref.index+1):ndrugs)){
    for( j in 1:nd.new[k]){
      OR[j,k] <- exp(s.beta.1[k]*new.dose[k,j]+  s.beta.2[k]*f.new.dose[k,j])
      odds.drug[j,k] <- OR[j,k]*exp(Z)
      p.drug[j,k] <- odds.drug[j,k]/(1+odds.drug[j,k])
    }
  }
  gamma[ref.index] <- 0
  for (k in c(1:(ref.index-1),(ref.index+1):ndrugs)){
    gamma[k] ~ dnorm(g,prec.g)

  }
  g ~dnorm(0,0.001)
  tau.g~dunif(0,5)      # common standard deviation is given a vague prior
  prec.g<-pow(tau.g,-2)
}












# mydoseresNMA <- function(){
#   b1[ref.index] <- 0
#   b2[ref.index] <- 0
#   for (i in 1:ns) {
#      rr[i,1]~dbinom(p0[i],nn[i,1])
#
#     # z[i] ~ dnorm(0,0.001)
#     # z[i] ~ dnorm(Z, prec.z)
#     delta[i, 1] <- 0
#     #d[i, 1] <- 0.00000E+00
#     u[i] ~ dnorm(0, 0.001)
#     for (k in 1:na[i]) {
#       r[i, k] ~ dbin(p[i, k], n[i, k])
#       logit(p[i, k]) <- theta[i,k]
#       theta[i,k] <- u[i] + delta[i, k]
#       rhat[i, k] <- p[i, k] * n[i, k]
#       dev[i, k] <- 2 * (r[i, k] * (log(r[i, k]) - log(rhat[i,k])) + (n[i,k] - r[i, k]) * (log(n[i, k] - r[i,k]) - log(n[i, k] - rhat[i, k])))
#     }
#     resdev[i] <- sum(dev[i, 1:na[i]])
#     for (k in 2:na[i]) {
#       d[i, k] <- (b1[t[i, k]] * dose1[i, k]) - (b1[t[i, 1]] * dose1[i, 1]) + (b2[t[i, k]] *dose2[i, k]) - (b2[t[i, 1]] * dose2[i,1])
#       delta[i, k] <- d[i, k]
#     }
#   }
#
#   for (k in c(1:(ref.index - 1), (ref.index + 1):ndrugs)) {
#     b1[k] ~ dnorm(0.00000E+00, 0.001)
#     b2[k] ~ dnorm(0.00000E+00, 0.001)
#   }
#
#
#   for (m in 1:ns){ ## for each study
#  #   # rr[m,1] ~ dbinom(p0[m],nn[m,1])
#     logit(p0[m]) <- OO[m]
#     OO[m] ~ dnorm(Z, prec.z)
#  #   # z[m] ~ dnorm(0,0.001)
#   }
#  # priors
#  Z ~ dnorm(0, 0.001)
#  prec.z <- pow(sigma.z,-2)
#  sigma.z ~ dunif(0,10)
#
#  for (k in c(1:(ref.index-1),(ref.index+1):ndrugs)){
#    for( j in 1:nd.new[k]){
#      OR[j,k] <- exp(b1[k]*new.dose[k,j]+  b2[k]*f.new.dose[k,j])
#      odds.drug[j,k] <- OR[j,k]*exp(Z)
#      p.drug[j,k] <- odds.drug[j,k]/(1+odds.drug[j,k])
#    }
#  }
#  totresdev <- sum(resdev[])
# }
htx-r/doseresNMA documentation built on Jan. 28, 2021, 5:32 a.m.